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Unconventional phases often occur where two competing mechanisms compensate. An excellent 
example is the ionic Hubbard model where the alternating local potential <5, favoring a band insulator 
(BI), competes with the local repulsion U, favoring a Mott insulator (MI). By continuous unitary 
transformations we derive effective models in which we study the softening of various excitons. The 
softening signals the instability towards new phases that we describe on the mean-field level. On 
increasing U from the BI in two dimensions, we find a bond-ordered phase breaking orientational 
symmetry due to a d-wave component. Then, antiferromagnetic order appears coexisting with the 
d-wave bond order. Finally, the d-wave order vanishes and a Neel-type MI persists. 

PACS numbers: 71.30.-|-h,71.10.Li,71.10.Fd,74.20.Fg 


Searching for unconventional states of matter and non¬ 
trivial elementary excitations (quasiparticles (QPs)) is 
one of the crucial objectives of the research in strongly 
correlated lattice models. Outstanding examples range 
from the quasi long-range ordered Mott insulator (MI) in 
one dimension (ID) where the neutral spin-1/2 particle 
“spinon” represents the elementray excitation to quan¬ 
tum spin ice in three dimensions (3D) where magnetic 
monopoles represent the QPs. In order to find unex¬ 
pected phases, it is a good idea to focus on parameter 
regions where two antagonists compensate because then 
subtle subleading mechanisms can take over. 

The present article addresses the IHM in two dimen¬ 
sions (2D) on the square lattice in order to understand 
which phases possibly arise between the band insulator 
(BI) and the MI. Its Hamiltonian reads 

r,(7 r 

+ t X (4.<TCr.'.<7 +h-C.), (1) 

{rr') ,<7 

where r := ix + jy spans the square lattice, cj. ^ ^) 

is the fermionic creation (annihilation) operator at site 
r with spin cr, and g. := c]) a-Cr.a- the occupation 
operator. The sum over {rr') restricts the hopping to 
nearest-neighbor sites. Initially, the IHM was introduced 
to describe the neutral-ionic transition in mixed-stack or¬ 
ganic compounds^ and was later used for the description 
of ferroelectric perovskites^. 

In ID, it is well understood that the BI at small Hub¬ 
bard interaction U is separated from quasi-long-range or¬ 
dered MI at large U by an intermediate phase with al¬ 
ternating bond order (BO)^. The position of the two 
transition points {Ud from BI to BO and Uc 2 from BO 
to MI)^i^ and the excitation spectrunt^i^ of the model are 
determined quantitatively. We highlight that the transi¬ 
tion to the BO phase is signaled by the softening of an 
excitoit^Ti^ located at momentum tt (setting the lattice 
constant to unity )^ii. 


In 2D, the limiting BI and MI phases are expected at 
low and high U, respectively. The MI at large U shows 
long-range antiferromagnetic order of Neel type (AF)^. 
But existence and properties of intermediate phases is 
highly controversial. A determinant quantum Monte 
Carlo (DQMC) analysis of small clusters identifies an 
intermediate metallic region without magnetic or bond 
order—. Cluster dynamic mean-field theory (DMFT) 
indicates an intermediate phase with an incommensu¬ 
rate bond order and a finite staggered magnetizationiS. 
Single-site DMFT of the half-filled IHM on the Bethe lat¬ 
tice suggests a single transition point separating the BI 
from the magnetically ordered MI— i. 
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Figure 1. (color online) Scheme of the phase transitions in 
the 2D IHM. The staggered magnetization is denoted by 
and the spin-dependent d-wave bond order parameter by d^. 
The S' = 0 exciton gap is Ae and As is the spin gap (or S = 1 
exciton)^. 


We study the 2D IHM by continuous unitary trans¬ 
formations (CUTs)i^ realized in real space up to higher 
orders in the hopping t. The flow equations are closed 
using the directly evaluated enhanced perturbative CUT 
(deepCUT)i^. In this way, an effective model is de¬ 
rived in terms of the elementary fermionic QPs. The 
virtual processes are eliminated as in the derivation of 
a t-J model from a Hubbard modeliAJi. The effective 
model allows us to address the interaction of two of these 
QPs. A particle and a hole attract each other and may 
form a bound state, an exciton. The softening of such 
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excitons indicate the instability of the phase towards a 
condensation of these excitons. In this fashion, we deter¬ 
mined the BI-to-BO transition value Uci of the ID IHM 
within 0.3%^. We stress that already the second order 
calculation in the hopping t yields qualitatively the cor¬ 
rect result, i.e., the correct symmetries for the condensate 
phase, in spite of the very small energy scales /U (rel¬ 
ative to 5) driving the transition^. Thus we adopt the 
same strategy in 2D where the CUT has been successfully 
used befor o^^d^ , in particular for binding phenomenai 2 £. 

We are using the double 1-particle gap Ac := — 

Eq), the singlet exciton gap Ag := E^g^g — E^, and the 
spin (triplet) exciton gap Ag := E^g^^ — Eq where Eq 
stands for the ground state energy and E^g for the first 
excited state with total spin S for N electrons. At half¬ 
filling, the number of electron N equals the number of lat¬ 
tice sites. In the absence of any electron-electron bound 
state, Ac, being twice the minimum of the fermion dis¬ 
persion, equals the charge gap. We use the term charge 
gap because it has been used in previous papers on the 
IHMi-i. 

Fig. □ summarizes our key results: Upon increasing 
U in the BI, the exciton gap Ac reduces till it vanishes 
at Uci where d-wave bond order (dBO) sets in. In this 
phase, the 90° rotation symmetry of the square lattice is 
broken and the d-wave bond order parameter do- becomes 
finite. Then, the spin gap Ag closes at Uc 2 beyond which 
the staggered Neel magnetization mg becomes finite and 
the d-wave order parameter da starts depending on spin 
(dAF). Eventually, the rotational symmetry is restored 
at Uc 3 {da = 0) while the Neel order persists (AF). 

Technically, we proceed as in ID^ by transforming 
the half-filled IHM into a translationally invariant prob¬ 
lem. On the odd sublattice {i + j odd), we apply the 
electron-hole transformation : cj. ^ ^. For 

consistency, we then use for particle and hole creation 
the same operator The local phase transformation 
r('): e^^{^+y)e-^T yields 
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Figure 2. (color online) (a) Low-energy excitations for t = 
0.055 and U = 1.15 in order 6 ; charge gap Ac (b) and singlet 
gap As (c) in units of 5 versus U. 
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disentangle the one-QP and the two-QP subspaces from 
the remaining Hilbert space by a suitable change of ba¬ 
sis such that the 1-QP dispersion and the dispersion of 
bound states such as excitons become accessible. Techni¬ 
cally, we employ the reduced generator n : 321 targeting 
the 0-, 1-, and 2-QP sector—. The sign of the terms in 
the generator, i.e., the sense of rotation, is chosen first ac¬ 
cording to the change in the number of double occupancy 
(increase of double occupancy leads to the same sign in 
the generator as in the Hamiltonian). If this number is 
not changed, the change in the number of fermions de¬ 
termines the sign. In this way, we map the Hamiltonian 
(I^ to the effective Hamiltonian 


iL=(C/-25)/4^1+t ^ {Aafl,a + ^-0-) 

r <rr'>,cT 

+ (5 - U)/2 ^ + U ^ fyijrJr^V (2) 

r.a r 

The phase transformation induces a shift in momen¬ 

tum space by 7r/2 := (7r/2,7r/2) for single /-fermions. 
The electron-hole transformation pU'*) modifies the spin 
and the charge on the odd sublattice according to 


{S-, sy, {-s^,+sy, l - (3) 


Hence, the AF phase in the original c-fermions ap¬ 
pears as translationally invariant ferromagnet in the /- 
fermions. The staggered magnetization reads rUg := 
1 sgn(o')(/o afo a)- The d-wave order parameter reads 

(7 


da ■■= {fLfU - (/L/L)- 


We apply the deepCUT to ([2]) eliminating processes 
which let single /-fermions and pairs of them decay. We 


= + + h.c.) + Hi (4) 

r,a d 

where stands for interaction terms of two creation 
and two annihilation operators with coefficients C|. The 
effective coefficients Cl and C| are determined from the 
flow equation d{,H = [ 77 , ff] 12 , 13 , Higher particle in¬ 
teractions are ignored. The range of processes in (|1]) is 
limited by the order in t of the deepCUT. The calcula¬ 
tions are restricted to order 6 because for higher orders 
overlapping continua prevent a well-defined convergence 
of the flow2i. Thus the effective Hamiltonian is not fully 
quantitative, but it is qualitatively correct. For an assess¬ 
ment of the accuracy see Refs.@ and[3; technical aspects 
related to simplification rules are given in Appendix 
The effective Hamiltonian allows us to determine 
the ground state energy, the 1-QP dispersion, and bound 
pairs of 2 QPs. Due to the particle-hole transforma¬ 
tion excitons with odd distance between their con- 
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Figure 3. Correlated hopping Re¬ 

call that the spin configuration shown is a particle-hole singlet 
in the original c-fermions. 


ping, see Fig. [31 Such terms, occurring already in the 
second order in t with positive sign, have a singlet on a 
bond in i-direction hop to a bond in y-direction^. This 
correlated hopping corresponds to the well-known three- 
site terms in the t-J modeli^^—. If the wave function has 
the opposite sign on y-bonds of those on x-bonds this cor¬ 
related hopping leads to attraction and eventually to a 
bound state with d-wave symmetry. 


stituents appear in the effective model as pairs of par¬ 
ticles. The eigenvalues in the 2-QP subspace are found 
for fixed total momentum k, total spin S', and total 
component M as in the ID case^. Fig. HJa) depicts the 
low-energy spectrum in the BI phase consisting of the 
fermionic dispersion and a singlet exciton; no triplet ex- 
citon is found. Note the rather flat dispersion and lower 
continuum edge. 

Very interestingly, the singlet exciton takes its mini¬ 
mum at fc = (0,0). In ID, its minimum was at fc = tt 
so that its condensation led to dimerized bond orderi^i^. 
The fact that the singlet exciton softens at a different 
momentum clearly shows that the condensed phase be¬ 
yond the BI will no longer be a dimerized BO. Figs.|3](b) 
and (c) display the charge and the singlet gap vs. U. 
Generally, the softening of a bosonic excitation implies 
the divergence of the susceptibility with the same sym¬ 
metries x(w) at w = 0 because x(^) corresponds to the 
bosonic propagator up to matrix elements. Thus it has 
a pole at the boson eigen energy moving to tu = 0 upon 
softening. The direct access to the energies of bound 
states avoids the necessity to guess the correct symmetry 
of the diverging susceptibility. In order 4, the singlet gap 
Ae vanishes at Ud = 1.08(5 where the charge gap Ac is 
still finite at 0.003(5. In order 6, we find Ud = 1.126i5 
with Ac = 0.023(5. Considering all available orders and 
comparing to the ID results^, we estimate that the tran¬ 
sition point Ud lies between U = 1.08(5 and U = 1.126(5. 
The value of the charge gap at the transition point is 
also expected to be between the results in order 4 and 
order 6, similar to ID findings^. Increasing t to 0.1(5 the 
transition point shifts to larger values of [/, see Appendix 

E 

The key observation is that the S = 0 exciton soft¬ 
ens before the system would become metallic in analogy 
to the ID case, but with other symmetries. The exciton 
condensation leads to a phase without magnetic order be¬ 
cause the condensing exciton has no spin. The condensa¬ 
tion does not break translational invariance because the 
exciton condenses at zero total moment. The analysis of 
the point group properties of the exciton reveals that its 
wave function becomes negative upon spatial rotation by 
90° so that the condensed phase should display the same 
d-wave symmetry. 

The mechanism for electron and hole to form a bound 
state with d-wave symmetry relies on the interaction 
fl+S:dfl+x+y,afr+xdfr,a representing correlated hop- 
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Figure 4. (color online) Staggered magnetization density ms 
and d-wave bond order parameter da (panel (a)) and charge 
gap A^, singlet gap Ae, and spin gap As (panels (b) and (c)) 
in units of 5 versus U for t = 0.055 at order 4. 

Similar to the ID case^, we describe this condensed 
phase by a self-consistent BCS-type MF theory applied 
to dH), see Fig. HKa) to the left of Uc 2 - This approach is 
exact for U < Ud where the quantum fluctuations are 
already considered by the CUT. In particular, the onset 
of order occurs precisely where the exciton energy van¬ 
ishes as required by consistency. For U > Ud the MF ap¬ 
proach is a good approximation systematically controlled 
by the distance U > Ud ■ The order of the calculations is 
4 in powers of t and t = 0.05(5. Its bond order da consists 
in the difference of the bonds in x and in y direction im¬ 
plying d-wave character. We determine the expectation 
values ifl^afdd') (/q^r/d.a) self-consistently, see Ap¬ 
pendix jC] comprising also the magnetization rris and the 
d-wave order parameter da- 

The singlet and the triplet excitons within the dBO 
are determined by analysing the 2-QP problem neglect¬ 
ing the off-diagonal terms linking 2-QP states to states 
of higher particle number, see Appendix |C] for the de¬ 
tails. This leads to the gaps displayed in Fig. SKb). The 
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exciton gap Ae quickly increases till the exciton ceases 
to exist because it merges with the 2-QP continuum. In 
parallel, a triplet exciton is formed and softens. The re¬ 
sults of order 6, given in Appendix [HI are qualitatively 
similar. Hence the system becomes unstable towards the 
condensation of a magnetic bosons, i.e., towards mag¬ 
netic ordering. The softening occurs at momentum tt in 
the original model so that this triplon condensation is the 
expected Neel ordering. 

We study the Neel ordering by allowing for spin de¬ 
pendent MF solutions of ([3]), i.e., mg becomes finite and 
df ^ di- In addition, the charge gap becomes also spin 
dependent A| ^ A;t. We find that the d-wave order 
does not vanish, but that it persists so that antiferro¬ 
magnetism and d-wave order coexist and a d-wave anti- 
ferromagnet (dAF) forms. This solution occurs at and 
beyond Uci « 1.083(5 as shown in Fig. [DJa) and in the 
zoom (c)i^ 

So far, we constructed the phase diagram following the 
spontaneous symmetry breaking indicated by soft exci- 
tons. In the dAF phase, the symmetry is strongly re¬ 
duced. For large t/, however, we expect Neel ordering 
without d-wave components because the Hubbard model 
maps to the Heisenberg model since the ionic alterna¬ 
tion is only a small perturbation tor U ^ 6. The MF 
solution retrieves this feature. The BO parameters d^ 
quickly decrease on increasing U. For U > 1.123(5, d^ is 
smaller than 10“^. We conclude that there is a Ucs be¬ 
yond which the spatial rotation is restored as a symmetry 
of the IHM: the dAF phase becomes a Neel ordered AF 
phase. The precise value of Ucs is difficult to pinpoint 
because d^{U) is consistent with an exponential decrease 
cx exp((7/(1/ — Ucs)) with Ucs ~ 1.245. This resembles 
a Kosterlitz-Thouless type of transition as also proposed 
in Ref. M But the numerics is too involved and the en¬ 
ergy scales too small to identify the decrease of df{U) 
unambiguously. 

In Fig. EKc), we also plot the spin-dependent charge 
gaps in the vicinity of Uc 2 - Beyond Uc 2 , A/ increases 
quickly, but Aj decreases and hits zero at a half-metallic 
point C/hm- Such half-metallic points occurred in previ¬ 
ous analyses alreadyii. The spin-up charge gap is located 
at momentum fc = 7r/2 up to C/hm, but for U > C/hm it 
moves to incommensurate momenta at the line kx = ky 
and remains zero. But we cannot exclude that the 
zero A| for U > C/hm is an artefact of our approach. 
Any small deviation from the particle-conserving effec¬ 
tive Hamiltonian would spoil this feature. Thus we 
also included ™ ®c) displaying a more generic 

behavior. 

Summarizing, we studied unconventional phases in the 
ionic Hubbard model in 2D between the known band in¬ 
sulator at weak interaction and the Neel ordered Mott 
insulator at strong interaction. We found bond order 
which spontaneously breaks the orientational symmetry: 
d-wave bond order (dBO). Even a region of coexistence of 
bond order and Neel order is identified: d-wave bond or¬ 
der and antiferromagnetism (dAF). Thereby, we provide 


clarification for the so far ambiguous evidence concern¬ 
ing the influence of strong interactions on band insula¬ 
tors. We stress that no weak-coupling instabilities indi¬ 
cate the existence of these phases. This underlines that 
these phases are inherent to strong coupling^^. Against 
this background, we think that the single-site DMFT— 
does not capture the spatial correlations indispensable 
for d-wave ordering and the DQMCS considers too small 
clusters at too large temperatures to capture the small 
energy scales driving the transitions found. The clus¬ 
ter DMFT— also advocated bond order, but combined 
with dimerization and incommensurability. Since the ro¬ 
tational symmetry was not allowed to be broken no d- 
wave ordering could be detected. 

We emphasize that the scenario discovered here is dif¬ 
ferent from the bond ordered phases with d-wave sym¬ 
metry discussed in the literature so far. These pre¬ 
vious scenarios imply either a kind of alternation^i^l 
or incommensurability22i2^ and they usually occur away 
from half-filling in contrast to our finding. The dBO 
phase presently advocated induces no dimerization and 
breaks only the rotational symmetry of the square lattice. 

The conclusions are based on effective models ex¬ 
pressed in elementary excitations. This approach was 
successfnl in ID. No ad hoc assumptions about broken 
symmetries are made, but instabilities are systematically 
deduced from softening bound states, namely singlet and 
triplet excitons, the latter being the spin excitations. The 
quantum numbers and symmetries of the excitonic wave 
function implies the symmetry of the phase to which the 
instable phase evolves. In this way, two continuous phase 
transitions at Uci and Uc 2 are identified. Based on sym¬ 
metry considerations, we expect the first transition to be 
in the Ising and the second transition to be in the (7(3) 
universality class. 

We consider the obtained scenario to be complete be¬ 
cause the d-wave order quickly decreases in the dAF 
phase so that only Neel order remains beyond a third 
transition interaction Ucs as expected for large interac¬ 
tion U. 

The findings indicate for which types of correlations 
one has to look in order to identify intermediate phases 
in the wide research field between band insulating and 
Mott insulating phases, for instance on various lattices. 


ACKNOWLEDGMENTS 

We are grateful to the Helmholtz Virtual Institute 
“New states of matter and their excitations” for finan¬ 
cial support. 


Appendix A: Simplification Rules 

Here we present the simplification rules used to realize 
the comprehensive CUTs, cf. Ref. [l^. 
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The basic a-posteriori and a-priori simplification rules 
(SRs) introduced in Ref. H can be used for the IHM 
on the square lattice because they are lattice indepen¬ 
dent. In addition, we have implemented an extended 
a-posteriori SR. The aim is to estimate whether a mono¬ 
mial A contributes to the targeted quantities up to a 
specific order, for example n, or not. We explain the ex¬ 
tended a-posteriori SR for the creation operators in A. 
The annihilation part can be analyzed in the same way. 

Similar to the ID case, one can treat spin-up and spin- 
down operators separately^. From the position of cre¬ 
ation operators with spin cr, we form a cluster (graph). 
We consider a vertex for each lattice site and add an 
edge between two vertices if the corresponding sites are 
adjacent. The degree of a vertex is defined as the num¬ 
ber of edges linked to the vertex. The cluster is divided 
into a set of linked subclusters. Let us consider one of 
the subclusters and denote it by C. In the first step, we 
need to underestimate the number of commutations K[C\ 
necessary to cancel the arbitrary linked cluster C. 

The size of C is reduced in an iterative way. We elimi¬ 
nate each vertex with degree 1 and the vertex it is linked 
to. This costs one commutation for each of these pairs 
to be omitted. Let us suppose that after commuta¬ 
tions no vertex with degree 1 is left. This reduces the 
linked cluster C to isolated vertices and/or to the linked 
subclusters C which has no vertex with degree 1. Each 
isolated vertex requires one commutation to be canceled. 
We denote by the number of isolated vertices. One 


also needs 


commutations to cancel the subclus¬ 


ter C which has v^., vertices. |"x] is defined as smallest 
integer not less than x. Hence, we obtain 


K[C\ > K [C] ■= p^ + Sc + ^ 


c 


> Pc + Sq + 


(Ala) 

(Alb) 


where Vc simpler Eq. (lAlbD 

c 

instead of (lAlaD for K [C]. There should be at least two 
subclusters C' each one with an odd number of vertices 
to change the “equal” sign in (lAlbl) into “greater than”. 
For the square lattice, this means that the initial linked 
cluster CF should contain at least 16 (=7-|-7-|-2) vertices. 

One can check that Eq. (|A1I) simplifies to the ID coun¬ 
terpart where no subcluster C' can exisli^. As example, 
let us examine the three clusters in Fig. [5] We find 
Pc = 1, Sc = 2, and no C for Fig. [DJa), Pc = 2, Sc = 2, 
and no C' for Fig. EKb), and finally Pc = I, Sc = 0, and 
one subcluster C' with vc = 4 for Fig. ISc). 

The total number of commutations Kq required to can¬ 
cel all the creation operators of A is given by 



(a) (b) (c) 

Figure 5. Three kinds of linked clusters on the square lattice. 


where the sum runs over all linked spin-up and spin-down 
clusters of monomial A. Similar to the rules presented 
in Ref. M Eq. dSH) can be extended to where sec¬ 
tors up to q quasiparticles are targeted. The aim is to 
keep q operators such that the maximal order of A is 
overestimatedi^. We first keep the operators which can 
save one commutation. The total number of these oper¬ 
ators is given by 

(®C + mod 2) , (A3) 

c c 

where we suppose that (lAlbl) holds for K[C]. The re¬ 
maining operators (if q > d) can save one commutation 
per pair. Therefore, we find 

Kg = max ^Kq — 

where d' := min((i, q). In the same way, one can analyze 
the annihilation part of the monomial A to find Kg . The 
monomial A can be discarded if 

0„,in(A) >n-iF“-iF“, (A5) 

where Omin(A) is the minimal order of A^. 



Appendix B: Results for Larger Hopping 

In order to show that the results obtained in the main 
text also apply for larger values of the hopping. Fig. 
(HKa) represents the same quantities as Fig. ^a) but for 
t = 0.1(5 and U = 1.25^. Figs. Elb) and (c) are the same 
as Figs. HKb) and (c), but for t = 0.15. One can see that 
by increasing t from t = 0.055 to t = 0.15, the position of 
the first transition point is shifted to the larger values of 
U. The value of the charge gap at the transition point is 
also increased compared to the case for t = 0.055. This 
behaviour is similar to what is found for the ID IHM 
by DMRG^. Larger hopping increases the energy scale 
on which the relevant physics takes place in the ionic 
Hubbard model (IHM). 


Appendix C: BCS-Type Mean-Field Analysis 


AS:=^A[C] 

c 


(A2) 


The Hamiltonian 0 is normal-ordered with respect 
to the bilinear part using Wick’s theorem. The bilin- 
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Figure 6. (a) The same as Fig. [2ja) in the main text, but for 
t = 0.1(5 and U = 1.255. (b) and (c) are also the same as Fig. 
lab) and (c) in the main text, but for t = 0.15. 
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ear Hamiltonian is diagonalized in momentum space by 
a Bogoliubov transformation. The relevant expectation 
values are given by 


i / ^cos(fe-d)dfc, (Cla) 
(/d.a/L) = / ^sin(fc-d)dfc, (Clb) 

where := eg + 2 ^ cos(fc-d), := 2 ^ Fg sin(fe-d), 

d d 

and UJ^ := ■ 

The coefficients and F^ are defined as prefactors of 
normal-ordered hopping and Bogoliubov operators over 
the distance d, respectively. Conservation of total charge 
guarantees that no hopping term over odd and no Bogoli¬ 
ubov term over even distances occurs. One can omit the 
spin indices from Eqs. (EB) because the SU(2) symme¬ 
try is preserved in the d-wave bond order (dBO) phase. 
The general spin-dependent form helps us to analyze also 
phases with broken spin symmetry. After some standard 
calculations, the Hamiltonian (3]) be rewritten as 


Hefl 


EqI + ^-f ^ ^ 

fe.cT kqiQf ’<^2 


/32/3ir^2]fcq/Xl 5 5 

(72(71 L^2j ^fc/24-q/ ,^2 '^kj2 — qf,0i '^kl2-\-qi ,cr2 ^kl2—qi,<7i 


(C2) 


where the operator is defined such that it diagonal¬ 
izes the non-interacting part of the Hamiltonian. The co¬ 
efficients of the interaction potential in momentum space 
C| are given in terms of the real space interactions C| 
and the hopping coefficients and the Bogoliubov coef¬ 
ficients Fjj,. The “• • • ” stand for the quartic off-diagonal 
interactions Hf oc A'^A'^A^A and Hq oc A^A^A'^A'^ and 
their hermitian conjugates. 

Taking the ground-state energy of the mean-field ap¬ 
proximation as the true ground-state energy neglects the 
terms Hg + H^. They are systematically controlled by 
U — Uci, i.e., they are indeed small in the parameter 
regime considered. 

Considering the 1-particle dispersion uJk in Eq. (IC2D 
and/or the interactions in the 2-particle sector above the 
mean-field solution, we also neglect the terms Hf -|- iJg 
which are again systematically controlled by U—Uci- The 
singlet and the triplet 2-particle bound states in Eig. I^b) 
of the main text are determined by exact diagonalization 
in the 2-particle subspace for fixed total momentum fc, 
total spin S, and total magnetic number M. The remain¬ 
ing quantum number is the relative momentum q which 
is not conserved. 

There is, however, a subtle point to be noticed. In 


the 2-particle state := |0)> 

states with two electrons (two holes) and one electron 
and one hole are mixed. In order to stay in the half-filled 
case, we define |fc;q)^=° := ^(1^; + |fc; 

which has the net total charge Q = 0. The relative 
momentum q is now restricted, e.g., to 0 < < tt and 

-^ < qy < We employ triplet states {S = 1) with spe¬ 
cific polarizations rather than specific magnetic numbers. 
The triplet states are given by \tx) = '^(It/) + l/t)): 

l^y) = - lit)), - lU)), and the 

singlet state by |s) = ^(|tt) + lU))- 

Constructing the Hamiltonian matrix in momentum 
space is time consuming compared to the real space di¬ 
agonalization. This is the case because for each and 
Qf in (IC2l) . one has to perform three summations over 
weighted real space interaction coefficients to calculate 
/32/3ir^21 ^1f 

0-2 (TiL 2jq^ 


Appendix D: Results in Order t® 

In order to underline that the analysis in the main 
text is supported by results in higher order of the hop- 
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ping. Figs. Elja) and (b) represent the same quantities REFERENCES 

as Fig. Hl^a) for U < Uc 2 and Fig. [Hb), but in order 
6. Divergence of the flow equations prevents us to con¬ 
sider interactions beyond U = 1.155(5, so that we can not 
access the full phase diagram. 



Figure 7. The same as in Fig. |4ja) for U < Uc 2 and Fig. Sib), 
but in order 6 of the hopping. 
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